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Tree  mortality  has  traditionally  been  evaluated  in  forest  inventories  through  summaries  of  dead  trees  by 
location,  species,  and  several  causal  agents.  Although  these  methods  were  most  commonly  used,  in  order 
to  assess  forests  sustainable  development,  they  have  had  limited  use  in  detecting  mortality  trends  and 
development  dynamics.  This  study  proposes  the  application  of  survival  analysis  for  the  purpose  of 
analyzing  tree  mortality.  Individual  tree  growth  increments  were  used  to  estimate  survival  and  hazard 
functions  for  the  Elatia  forest  (Drama,  Northeast  Greece).  These  estimates  provided  indications  of 
regional  mortality  by  diameter  at  breast  height  (DBH)  and  diameter  growth  ( ADBH)  between  successive 
measurements.  Comparisons  of  survival/hazard  curves  and  tests  of  effects  of  species  and  crown  class  (CC) 
on  individual  survival  curves  were  conducted.  Survival  analysis  technique,  by  using  the  variables  of  DBH 
and  ADBH,  could  help  foresters  to  evaluate  regional  tree  mortality  trends,  and,  consequently,  forests 
sustainable  development. 
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1.  Introduction 

Tree  mortality  in  forest  inventories  has  traditionally  been 
evaluated  using  relatively  simple  summary  statistics.  Mortality 
information  available  to  foresters  has  typically  included  losses  in 
timber  volume  due  to  mortality,  summaries  of  mortality  causal 
agents,  spatial  locations  of  mortality,  and  mortality  trends  by  species 
[38].  Currently,  most  forest  management  reports,  in  order  to  assess 
sustainable  development,  use  this  mortality  evaluation  methodolo¬ 
gy  [1,2].  More  in-depth  mortality  analysis  has  historically  only  been 
facilitated  through  development  of  individual  tree  mortality  logistic 
mdels,  a  technique  that  has  had  limited  use  in  national  inventories 
[3,4]  and  may  be  inadequate  for  broadly  defining  forest  develop¬ 
ment  dynamics  [5].  Although  remotely  sensed  information  and 
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geographic  information  systems  have  greatly  aided  analysis  of  forest 
mortality,  the  basic  analytical  methodologies  of  forest  mortality 
evaluation  have  only  slowly  evolved  [6]. 

The  forest  sciences  have  historically  focused  on  developing 
individual  tree  mortality  sub-models  for  incorporation  into  growth 
and  yield  models  [7-11].  Although  other  sciences  that  monitor 
populations  of  living  organisms,  such  as  the  veterinary  and  medical 
sciences,  have  developed  methodologies  to  evaluate  mortality 
beyond  that  of  the  individual,  the  forest  sciences  have  relatively 
few  methodologies  for  evaluating  tree  population  mortality  and 
sustainable  development  assessment  [6].  Commonly  used  forest 
mortality  analytical  techniques  lack  methodology  for  incorporat¬ 
ing  the  time-dependent  nature  of  tree  mortality,  hypothesis 
testing,  censoring  of  observations,  and  tests  for  effects  of 
covariates.  Given  the  diseases  and  epidemics  that  have  greatly 
altered  forest  ecosystems  and  future  forest  health  issues  that  may 
occur,  techniques  for  evaluating  tree  mortality  and  forest  decline 
would  benefit  forest  scientists  and  managers  alike. 
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Analytical  methods  developed  by  the  medical  sciences, 
collectively  termed  survival  analysis,  may  provide  the  basis  for 
development  of  new  forest  mortality  analytical  techniques  [12,13]. 
Survival  analysis  is  most  often  defined  as  a  class  of  statistical 
methods  for  studying  the  occurrence  and  timing  of  events,  most 
often  death  [14,15,40].  Survival  analysis  is  unique  in  that  it  allows 
for  censoring  of  observations  (lack  of  exact  time  of  death)  and 
inclusion  of  time-dependent  covariates,  in  addition  to  dealing  with 
non-normal  distributions  [36,40].  Waters  [16]  first  proposed  using 
survival  analysis  to  address  forest  mortality  issues,  but  such 
applications  have  been  restricted  mainly  to  forest  inventories  in 
even-aged  forest  plantations  [17-20],  forest  research  plots  [21-23] 
and  stand  table  projections  [24].  Although  there  has  been 
significant  work  on  the  development  of  estimation  procedures 
for  survival  and  hazard  functions  in  forest  research  plots  [23,19] 
and  in  relation  to  diameter  growth  [25,26],  survival  analysis 
techniques  have  not  been  widely  applied  to  forest  inventory 
analyses  due  to  the  inherent  lack  of  detailed  time  and  age 
information  for  forest  inventories  [27].  Given  the  current  lack  of 
baseline  forest  inventory  mortality  analyses  techniques  [28]  and 
the  potential  that  survival  analysis  offers  [12,6,13],  an  examination 
of  survival  analysis  in  the  context  of  a  forest  inventory  is  warranted 
and  may  refine  analysis  of  tree  mortality  and  forests  sustainable 
development. 

The  primary  goal  of  this  study  is  to  estimate  and  interpret  the 
central  functions  of  survival  analysis  (Kaplan-Meier  survival  and 
hazard  functions),  on  a  time  scale  defined  by  growth  in  diameter  at 
breast  height,  for  Elatia  forest  (Drama,  Northeast  Greece). 
Specifically,  the  objectives  of  this  study  are: 

1.  to  use  diameter  at  breast  height  (DBH)  and  growth  in  diameter 
at  breast  height  (ADBH)  in  applying  survival  analyses  techni¬ 
ques; 

2.  to  determine  if  Kaplan-Meier  survival/hazard  functions  can 
represent  actual  mortality  trends  in  a  manner  practical  for 
ecological  interpretation; 

3.  to  determine  species  and  social  status  effects  to  mortality 
trends. 


2.  Materials  and  methods 

2.2.  Data 

The  data  for  this  study  came  from  two  successive  measure¬ 
ments  at  Elatia  forest  (Drama,  Northeast  Greece),  in  the  context  of 
students  practice  (Dimokritos  University,  Greece,  Department  of 
Forestry  and  Management  of  the  Environment  and  Natural 
Resources).  Sample  trees  were  surveyed  in  2005  and  re-measured 
in  2009  (Table  1).  Individual  trees  (observations)  were  included 
that  met  the  following  criteria:  alive  at  time  one  (2005)  and  either 
dead  or  alive  at  time  two  (2009),  and  DBH  >12  cm  at  time  one. 
Individual  tree  attributes  measured  at  time  one  that  were  included 
as  predictors  of  mortality  in  this  study  were  diameter  at  breast 
height  (DBH,  cm)  and  crown  class  (CC).  If  a  tree  was  dead  at  time 
two  then  its  DBH  was  set  equal  to  the  DBH  at  time  two  or  the  DBH 
at  time  one,  whichever  was  larger.  Since  a  tree’s  DBH  may  shrink 
following  death,  an  estimate  of  the  maximum  DBH  the  tree 
attained  before  death  would  better  benefit  survival  analysis  than 


Table  1 

Mortality  of  sample  trees  at  Elatia,  2005-2009. 


Species 

Number  of  measured  trees 

Number  of  trees  that  died 

Picea  abies 

142 

18 

Fagus  sylvatica 

1 

0 

Pinus  sylvestris 

11 

3 

an  estimate  of  a  decaying  bole  diameter.  CC  is  a  measure  of  a  tree’s 
dominance  in  relation  to  adjacent  trees  in  the  same  stand  and  is 
coded  as  follows:  1,  open  grown;  2,  dominant;  3,  codominant;  4, 
intermediate;  5,  suppressed.  ADBH  was  calculated  as  the 
difference  in  DBH  between  time  one  and  time  two.  For  trees  that 
died  during  the  re-measurement  interval  (censored),  ADBH 
represents  less  of  an  average  growth  rate. 


2.2.  Survival  analysis 


Kaplan-Meier  survival  and  hazard  functions,  are  used  to 
quantify  the  probability  distribution  of  mortality  in  a  population 
[29].  The  survival  function  is  defined  as  [14,15,36]: 

S(t)  =  P(T  >  t)  (1) 

where  S(t)  is  the  probability  that  a  death  occurs  at  some  time  T  at 
least  as  great  as  time  t,  but  is  not  constrained  except  for  being 
greater  than  0. 

The  hazard  function  is  an  instantaneous  mortality  rate  and 
hence  is  a  conditional  probability  defined  as  [15,36]: 


lim 

At-»0 


P(t  <T<  t  + At|T>t) 
At 


where  h(t)  is  the  probability  that  death  occurs  exactly  at  time  t, 
given  that  it  has  not  occurred  before  then. 

The  survival  function  may  be  estimated  non-parametrically  by 
using  the  life-table  method  given  by  [15]: 


m = fjo  -  hj)  (3) 

j= i 

where  for  interval  i,  q  is  the  start  of  time  and  hi  is  the  conditional 
probability  of  death.  For  i  =  1  and  hence  q  =  0,  the  survival 
probability  is  set  to  1.0. 

The  life-table  non-parametric  estimate  of  the  hazard  function  at 
the  midpoint  of  each  time  interval  is  given  by  [15]: 

d 

h(tim)  =  b,(m  -  (Wj/2)  -  (di/2))  (4) 

where  for  the  ith  interval,  tim  is  the  midpoint,  dt  the  number  of 
deaths,  the  width  of  the  interval,  n,  the  number  of  individuals  at 
the  beginning  of  the  interval,  and  w,  is  the  number  of  cases 
censored  (exact  time  of  death  cannot  be  ascertained)  within  the 
interval.  Note  that,  the  survival  and  hazard  functions  are 
mathematical  functions  of  each  other;  given  one,  we  can  compute 
the  other. 

The  null  hypothesis,  that  the  survival  functions  are  the  same  for 
two  groups  of  individuals,  may  be  tested  by  the  non-parametric 
logrank  test  statistic  given  by  [30]: 

UL  =  £(d„  -  e,j)  (5) 

j=  1 

where  UL  is  the  summation  over  all  unique  event  times  (in  both 
groups)  and  there  are  a  total  of  r  such  times,  dy  is  the  number  of 
deaths  that  occur  in  group  1  at  time  j  and  ey  is  the  expected 
number  of  events  in  group  1  at  time  j.  The  expected  number  of 
events  is  given  by  nydjlnj ,  where  rq  is  the  total  number  of  cases  that 
are  at  risk  just  prior  to  time  j,  ny  the  number  at  risk  just  prior  to 
time  j  in  group  1,  and  dj  is  the  total  number  of  deaths  at  time  j  in 
both  groups.  Squaring  and  dividing  UL  by  the  estimated  variance 
provides  a  j2  statistic.  Additionally,  logrank  tests  may  be 
generalized  to  test  whether  covariates  are  associated  with  survival 
times. 

As  evidenced  in  survival  analysis  formulations,  time  to  an  event 
is  the  defining  component  of  survival  methods.  Hence,  the  major 
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obstacle  cited  as  limiting  the  application  of  survival  analysis  to 
forest  inventories  is  the  lack  of  specific  tree  ages  and  the  censoring 
of  tree  mortality  [27].  However,  knowledge  of  age  is  not  necessary 
for  implementation  of  survival  analyses.  Any  measurement  unit 
that  indicates  changes  in  an  individual’s  status  between  re¬ 
measurements  may  replace  the  traditional  survival  analysis 
variables  of  age  and  time.  For  forest  inventories  that  re-measure 
trees  at  regular  intervals,  e.g.,  national  inventories  [31,4],  DBH  and 
ADBH  may  assign  individual  trees  within  a  population  to  classes 
defined  by  tree  size  and  growth.  Whereas  medical  studies  may 
determine  survival  functions  for  demographic  cohorts  across 
calendar  years,  forest  inventory  survival  functions  may  be 
determined  for  DBH  classes  across  growth  [12]. 

In  this  study,  time  starts  at  the  first  forest  inventory,  when  a 
subject  begins  to  be  at  risk  for  the  event  or  begins  to  be  monitored 
for  the  event.  Stating  this  in  terms  of  DBH,  time  is  ADBH  (the 
increase  in  DBH  from  initial  survey).  Other  studies  [25,26]  have 
found  some  success  in  using  predicted  ADBH  to  refine  estimates  of 
time  of  tree  death.  Survival  function  in  this  study  S(ADBH)  gives 
the  probability  that  a  tree  will  continue  to  live  until  its  diameter 
has  increased  by  at  least  ADBH.  Hence,  for  all  previously  stated 
formulations  of  the  non-parametric  survival  and  hazard  function 
estimators  and  logrank  tests,  we  substituted  ADBH  for  time. 

Several  software  packages  produce  estimates  of  the  survival 
and  hazard  functions.  In  this  study,  the  SPSS  SURVIVAL  [30]  was 
used.  Trees  were  grouped  by  initial  DBH  into  10-cm  diameter 
classes.  The  survival  and  hazard  functions  (Eqs.  (1)  and  (2))  were 
compared  for  their  forest  science  applicability.  In  addition,  the 
survival  function  was  further  examined  in  terms  of  effects  of 
covariates  (DBH  class  and  CC),  using  logrank  tests  (Eq.  (5)). 

3.  Results 

The  survival  and  hazard  functions  were  estimated  separately 
for  the  three  initial  DBH  classes  for  the  three  species  in  the  2005- 
2009  period  (Figs.  1  and  2).  We  could  say  that  trees  of  the  middle 
diameter  class  ((22,32)  cm)  seem  to  have  the  highest  survival 
probability  (Fig.  1 ).  For  all  three  initial  diameter  classes,  the  hazard 
of  death  increases  while  ADBH  increases.  We  could  say  that  the 
smallest  trees  with  the  most  diameter  growth  had  the  greatest 
mortality  hazard  (Fig.  2). 

Hazard  functions  for  the  [  1 2,22  ]  cm  DBH  class,  which  appears  to 
have  the  greatest  mortality,  were  stratified  by  species  group,  in 
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Fig.  1.  Kaplan-Meier  survival  functions  for  time  one  diameter  classes  by  ADBH. 
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Fig.  2.  Hazard  functions  for  time  one  diameter  classes  by  ADBH. 
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Fig.  3.  Hazard  functions  for  time  one  DBH  class  [12,22]  cm  for  the  two  tree  species 
by  ADBH. 


Table  2 

Logrank  tests  for  homogeneity  of  survival  distributions  by  DBH  class  and  CC. 


Covariate 

x2 

Degrees  of  freedom 

Ti|jifj  p 

DBH  class  (cm) 

4.099 

2 

0.129 

Crown  class  (CC) 

0.915 

4 

0.922 

order  to  evaluate  specific  mortality  trends.  Risk  of  mortality  was 
distinctly  different  between  two  species  groups  ( Picea  abies  and 
Pinus  sylvestris),  across  the  biggest  values  of  ADBH  (Fig.  3). 

Logrank  tests  (Table  2),  which  were  conducted  for  effects  of 
covariates  on  survival  functions  (DBH  class  and  CC)  did  not  show 
any  significant  effect  (p  >  0.05). 

4.  Discussion 

A  longitudinal  unit  can  be  any  unit  that  measures  a  variable’s 
transition  from  one  state  to  another  [36].  The  greatest  hurdle  in 
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applying  survival  analytical  techniques  to  forest  inventories  is 
finding  appropriate  longitudinal  units  to  quantify  survival 
probabilities  and  mortality  hazards  [13].  If  time  or  ages  are  used 
as  longitudinal  units  in  forest  inventory  analyses,  a  number  of 
problems  may  be  encountered  as  evidenced  by  previous  work  in 
forest  survival  analysis  [12,22,23,19].  First,  all  observations  are 
censored.  The  exact  time  of  tree  death  is  uncertain  with  the 
inventory  re-measurement  date  often  serving  as  the  longitudinal 
measure.  Second,  the  survival  function  curve  is  partially  dependent 
on  when  and  where  the  measurements  were  taken.  For  example,  if 
the  bulk  of  mortality  is  located  in  a  certain  area  of  the  state  that  is 
inventoried  at  a  discrete  point  in  time,  then  the  resulting  survival 
curve  will  be  biased  if  time  is  used.  Third,  the  age  of  a  tree  is 
difficult  to  estimate  in  forest  inventories.  Previous  studies  in  even- 
aged  stand  conditions  have  been  able  to  develop  rather  flexible 
survival/hazard  functions  using  individual  tree  age  and  time 
[22,23,19].  However,  DBH  is  a  quantity  that  hypothetically 
increases  with  time  until  a  tree  dies  (Harcombe  [12]).  Although 
DBH  is  not  a  perfect  predictor  of  tree  age,  especially  in  uneven- 
aged  stands,  tree  diameter  could  be  used  as  a  hypothetical 
surrogate  for  age,  considering  that  ADBH  is  zero  at  time  one  and 
tree  survey  stops  at  time  two,  where  ADBH  >  0.  Time  (years)  may 
greatly  relate  to  the  survivalship  of  humans,  while  tree  growth 
over  intervals  of  time  (i.e.,  annual  diameter  growth)  may  be  a  more 
meaningful  metric  in  forest  ecology.  Bigler  and  Bugmann  [25] 
found  that  using  ADBH  models  to  refine  estimates  of  time  of  death 
greatly  improved  logistic  mortality  models.  Unfortunately,  the 
longer  the  inventory  re-measurement  interval  the  greater  the 
possibility  that  ADBH  may  inaccurately  reflect  time  due  to  the 
censoring  of  tree  death  and  variable  tree  growth  rates  (particularly 
in  uneven-aged  stands). 

There  are  numerous  estimation  procedures  for  the  survival  and 
hazard  functions,  each  with  its  own  advantages  and  drawbacks. 
This  study  used  the  elementary  life-table  Kaplan-Meier  estimation 
procedure  as  a  first  attempt  to  apply  survival  analysis  techniques 
to  a  forest  inventory;  a  technique  suggested  by  Harcombe  [12]  and 
Woodall  et  al.  [32].  The  Kaplan-Meier  approach  allowed  non- 
parametric  estimation  of  survivorship,  the  ability  to  compare 
survivorship  among  stratified  sample  units,  and  the  ability  to  test 
association  between  covariates  and  survivorship.  The  survival 
function  quantifies  mortality  cumulatively  through  the  diameter 
distribution,  while  the  hazard  function  may  display  specific  DBH 
midpoint  mortality  rates.  As  suggested  by  Manion  and  Griffin  [28], 
the  quantification  of  rates  of  mortality  across  diameter  classes 
helps  identify  atypical  mortality  trends  as  soon  as  they  arise.  The 
hazard  and  survival  functions  can  together  provide  an  initial 
evaluation  of  tree  mortality  and  assessment  of  sustainable 
development  for  forest  inventories  as  long  as  the  survey  interval 
of  time  is  approximately  the  same  between  re-measurements. 
Large  sampling  intervals  will  affect  ADBH  and  ultimately  the 
interpretation  or  application  of  syrvival  analyses. 

Traditionally,  insects  and  disease  tree  mortality  has  been 
expressed  in  terms  of  ratios  of  tree  mortality.  Hazard  functions 
can  be  used  for  more  detailed  analysis  of  mortality  dynamics  for 
any  tree  population  of  interest.  He  and  Alfaro  [33]  found  that 
survival  analysis  was  useful  in  analyzing  tree  resistance  to  pest 
attack;  however,  they  also  found  that  tree  survival  time  was 
related  to  seasonal  temperatures  and  precipitation.  Hazard 
functions  also  allow  a  broad  comparison  of  mortality  risk  rates 
among  species  and  diameter  classes.  For  forest  inventories, 
hazard  functions  may  aid  investigations  between  forest  mortality 
and  causal  agents  [12]. 

With  logrank  test  we  can  test  the  hypothesis  concerning 
survival/hazard  curves  comparison.  In  this  study,  the  logrank  test 
for  effects  of  covariates  on  the  survival  function,  time  one  DBH 
class  and  CC,  indicated  that  these  covariates  were  not  significant 


in  tree  survival.  As  found  in  other  mortality  studies,  crown 
conditions  may  be  an  important  predictor  of  individual  tree 
mortality  [34].  The  ability  to  associate  individual  tree  traits  with 
mortality  hazard  has  enormous  potential  benefits  for  the  study  of 
natural  processes  [13]. 

5.  Conclusion 

Forest  inventory  mortality  analyses  have  predominantly  been 
focused  on  logistic  regression  modelling  at  the  individual  tree- 
scale  and  simple  data  summarizations.  This  study  proposes  an 
approach  to  forest  mortality  evaluation  and  sustainable  devel¬ 
opment  assessment  involving  combination  of  established 
survival  modelling  techniques  (survival/hazard  functions)  with 
traditional  quantifications  of  forest  stand  attributes  (DBH 
distribution  and  diameter  growth).  There  was  an  attempt  to 
use  DBH  and  ADBH  as  surrogates  for  age  and  time  in  application 
of  survival  analysis.  We  found  that  there  are  no  significant 
differences  in  survival,  between  diameter  and  crown  classes,  for 
the  Elatia  forest.  P.  abies  seems  to  have  the  greatest  risk  of 
mortality  (compared  to  P.  sylvestris),  at  small  diameter  and  big 
diameter  growth.  This  study’s  approache  may  eventually  lead  to 
more  efficient  and  statistically  defensible  evaluation  of  tree 
mortality  and  sustainable  development  assessment  for  tree 
populations  across  different  forest  types,  locations,  and  suffering 
from  varying  damage  agents. 
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